Introduction

Performs a Bayesian calculation with uninformed priors to estimate model 2 of

\[ Y_{icb2} = \alpha_{b} + \delta^{T} T_{c} + \delta^{GK} T^{GK}_c + \beta X_{icb1} + \rho Y_{icb1} +\gamma_{1} \tau_{c} +\epsilon_{icb2} \]

Warnings (if applicable)

Setup

the needful

Load the relevant packages

Set global core settings, if applicable. Use global setting only if necessary, best to use as an inline setting to avoid over allocating system resources

The models - table 3

These models build out table 3 which presents the results for the primary and secondary outcomes utilizing informed priors with 1. normal distributions and 2. Cauchy distributions.We follow the advice of Rachael Meager and incorporate disagreement in the literature and use 6 standard deviations and a mean = 0

Per Capita Monthly Consumption - Primary Outcome

Load the Monthly per capita consumption model variables

Normal Distribution

per capita consumption: This is the basic benchmarking model utilizing the a global normal distribution with mean = 0 and 6 standard deviations.

Model Summery

summary(per_cap_consumption_normal_bayesmodel)
##  Family: gaussian 
##   Links: mu = identity; sigma = identity 
## Formula: consumption_asinh | weights(samp_wgt) ~ cost_deviation + treat_any + treat_GK + consumption_asinh_R1 + Lhh_wealth_asinh + Lvill_eligible_ratio + (1 | block) + (1 | vid) 
##    Data: per_cap_consumption_data (Number of observations: 1750) 
##   Draws: 4 chains, each with iter = 2000; warmup = 1000; thin = 1;
##          total post-warmup draws = 4000
## 
## Group-Level Effects: 
## ~block (Number of levels: 22) 
##               Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sd(Intercept)     0.27      0.06     0.16     0.40 1.00      925     1840
## 
## ~vid (Number of levels: 248) 
##               Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sd(Intercept)     0.38      0.03     0.32     0.44 1.00     1822     2416
## 
## Population-Level Effects: 
##                      Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS
## Intercept                8.51      0.18     8.16     8.87 1.00     4268
## cost_deviation           0.00      0.00     0.00     0.00 1.00     2799
## treat_any                0.11      0.08    -0.05     0.28 1.00     2023
## treat_GK                -0.13      0.08    -0.30     0.03 1.00     2313
## consumption_asinh_R1     0.18      0.01     0.16     0.21 1.00     6550
## Lhh_wealth_asinh         0.02      0.01     0.01     0.03 1.00     8404
## Lvill_eligible_ratio     0.16      0.33    -0.48     0.80 1.00     1889
##                      Tail_ESS
## Intercept                3113
## cost_deviation           2941
## treat_any                2617
## treat_GK                 3076
## consumption_asinh_R1     3162
## Lhh_wealth_asinh         2465
## Lvill_eligible_ratio     2542
## 
## Family Specific Parameters: 
##       Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sigma     1.16      0.01     1.14     1.19 1.00     6660     2748
## 
## Draws were sampled using sampling(NUTS). For each parameter, Bulk_ESS
## and Tail_ESS are effective sample size measures, and Rhat is the potential
## scale reduction factor on split chains (at convergence, Rhat = 1).

Prior summery - how informative are priors

prior_summary(per_cap_consumption_normal_bayesmodel)
##                    prior     class                 coef group resp dpar nlpar
##              normal(0,6)         b                                           
##              normal(0,6)         b consumption_asinh_R1                      
##              normal(0,6)         b       cost_deviation                      
##              normal(0,6)         b     Lhh_wealth_asinh                      
##              normal(0,6)         b Lvill_eligible_ratio                      
##              normal(0,6)         b            treat_any                      
##              normal(0,6)         b             treat_GK                      
##  student_t(3, 10.7, 2.5) Intercept                                           
##     student_t(3, 0, 2.5)        sd                                           
##     student_t(3, 0, 2.5)        sd                      block                
##     student_t(3, 0, 2.5)        sd            Intercept block                
##     student_t(3, 0, 2.5)        sd                        vid                
##     student_t(3, 0, 2.5)        sd            Intercept   vid                
##     student_t(3, 0, 2.5)     sigma                                           
##  bound       source
##                user
##        (vectorized)
##        (vectorized)
##        (vectorized)
##        (vectorized)
##        (vectorized)
##        (vectorized)
##             default
##             default
##        (vectorized)
##        (vectorized)
##        (vectorized)
##        (vectorized)
##             default
check_prior(per_cap_consumption_normal_bayesmodel)
##                Parameter Prior_Quality
## 1            b_Intercept uninformative
## 2       b_cost_deviation uninformative
## 3            b_treat_any uninformative
## 4             b_treat_GK uninformative
## 5 b_consumption_asinh_R1 uninformative
## 6     b_Lhh_wealth_asinh uninformative
## 7 b_Lvill_eligible_ratio uninformative

Diagnostics

# trace diagnostic plot
mcmc_trace(per_cap_consumption_normal_bayesmodel, n_warmup = 0,
           pars = c("b_Intercept", "b_cost_deviation", "b_treat_any", 
                    "b_treat_GK", "b_consumption_asinh_R1", "b_Lhh_wealth_asinh",
                    "b_Lvill_eligible_ratio", "sd_block__Intercept", 
                    "sd_vid__Intercept", "sigma"))

ggsave("table_3_diagnostics\\per_cap_consumption_norm_trace.png", plot = last_plot(), width = 12, height = 5)


#density diagnostic plot
mcmc_dens(per_cap_consumption_normal_bayesmodel,
           pars = c("b_Intercept", "b_cost_deviation", "b_treat_any", 
                    "b_treat_GK", "b_consumption_asinh_R1", "b_Lhh_wealth_asinh",
                    "b_Lvill_eligible_ratio", "sd_block__Intercept", 
                    "sd_vid__Intercept", "sigma"))

ggsave("table_3_diagnostics\\per_cap_consumption_norm_dens.png", plot = last_plot(), width = 12, height = 5)

mcmc_dens_overlay(per_cap_consumption_normal_bayesmodel,
           pars = c("b_Intercept", "b_cost_deviation", "b_treat_any", 
                    "b_treat_GK", "b_consumption_asinh_R1", "b_Lhh_wealth_asinh",
                    "b_Lvill_eligible_ratio", "sd_block__Intercept", 
                    "sd_vid__Intercept", "sigma"))

ggsave("table_3_diagnostics\\per_cap_consumption_norm_overlay.png", plot = last_plot(), width = 12, height = 5)

mcmc_dens(per_cap_consumption_normal_bayesmodel,pars = c("b_treat_GK"))

mcmc_dens_overlay(per_cap_consumption_normal_bayesmodel,pars = c("b_treat_GK"))

#acf (auto-correlation) diagnostic plot
mcmc_acf(per_cap_consumption_normal_bayesmodel,
           pars = c("b_Intercept", "b_cost_deviation", "b_treat_any", 
                    "b_treat_GK", "b_consumption_asinh_R1", "b_Lhh_wealth_asinh",
                    "b_Lvill_eligible_ratio", "sd_block__Intercept", 
                    "sd_vid__Intercept", "sigma"))

ggsave("table_3_diagnostics\\per_cap_consumption_norm_acf.png", plot = last_plot(), width = 12, height = 5)

posterior predictive checks

pp_check(per_cap_consumption_normal_bayesmodel, nsamples = 100)
## Warning: Argument 'nsamples' is deprecated. Please use argument 'ndraws'
## instead.

pp_check(per_cap_consumption_normal_bayesmodel, nsamples = 10, type = 'error_scatter_avg', alpha = .1)
## Warning: Argument 'nsamples' is deprecated. Please use argument 'ndraws'
## instead.

Dietary Diversity - Primary Outcome

Load the Dietary Diversity variables

Dietary Diversity: This is the basic bechmarking model utilizing the default, uninformed priors

dietary_diversity_normal_bayesmodel <- 
  brm(formula = dietarydiversity | weights(samp_wgt) ~ 
        cost_deviation + treat_any + treat_GK + 
        dietarydiversity_R1 + Lhh_wealth_asinh + Lvill_eligible_ratio + Lsavingsstock_asinh3 + 
        Lconsumpti_x_Ldietarydi + Lconsumpti_x_Lproductiv + Ldietarydi_x_Lassetscon + 
        (1 | vid) + (1 | block),
     data = dietary_diversity_data,
     family = gaussian("identity"),
     set_prior("normal(0,6)", class = "b"),
     seed = 1272022,
     warmup = 1000,
     iter = 2000,
     thin = 1,
     control = list(adapt_delta = .95, max_treedepth = 10),
     #backend = "cmdstanr",
     cores = 4, #overrides default 1 core
     #threads = 3,need to get cmdstanr package working here
     save_pars = save_pars(all = TRUE), # potentially allows for more post-processing functionality
     file = "informed_prior_outcomes\\dietary_diversity_normal_bayes")
# tidy_dietary_diversity_normal_bayesmodel <- tidy(dietary_diversity_normal_bayesmodel)
# #view(tidy_dietary_diversity_normal_bayesmodel)
# write_csv(tidy_dietary_diversity_normal_bayesmodel, "informed_prior_outcomes\\dietary_diversity_normal_bayes.csv")

Model Summary

summary(dietary_diversity_normal_bayesmodel)
##  Family: gaussian 
##   Links: mu = identity; sigma = identity 
## Formula: dietarydiversity | weights(samp_wgt) ~ cost_deviation + treat_any + treat_GK + dietarydiversity_R1 + Lhh_wealth_asinh + Lvill_eligible_ratio + Lsavingsstock_asinh3 + Lconsumpti_x_Ldietarydi + Lconsumpti_x_Lproductiv + Ldietarydi_x_Lassetscon + (1 | vid) + (1 | block) 
##    Data: dietary_diversity_data (Number of observations: 1751) 
##   Draws: 4 chains, each with iter = 2000; warmup = 1000; thin = 1;
##          total post-warmup draws = 4000
## 
## Group-Level Effects: 
## ~block (Number of levels: 22) 
##               Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sd(Intercept)     0.30      0.08     0.16     0.48 1.01      999     1490
## 
## ~vid (Number of levels: 248) 
##               Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sd(Intercept)     0.55      0.04     0.47     0.64 1.00     1524     1977
## 
## Population-Level Effects: 
##                         Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS
## Intercept                   3.19      0.22     2.75     3.62 1.00     3506
## cost_deviation              0.00      0.00     0.00     0.00 1.00     2913
## treat_any                   0.25      0.12     0.02     0.48 1.00     2215
## treat_GK                   -0.03      0.11    -0.26     0.19 1.00     1985
## dietarydiversity_R1        -0.12      0.07    -0.27     0.02 1.00     3123
## Lhh_wealth_asinh            0.01      0.01    -0.00     0.03 1.00     5196
## Lvill_eligible_ratio       -0.57      0.44    -1.46     0.30 1.00     1857
## Lsavingsstock_asinh3        0.00      0.00     0.00     0.00 1.00     4079
## Lconsumpti_x_Ldietarydi     0.01      0.01     0.00     0.03 1.00     3153
## Lconsumpti_x_Lproductiv     0.01      0.00     0.00     0.01 1.00     4034
## Ldietarydi_x_Lassetscon     0.01      0.00     0.01     0.01 1.00     4947
##                         Tail_ESS
## Intercept                   2970
## cost_deviation              2827
## treat_any                   2894
## treat_GK                    2504
## dietarydiversity_R1         2880
## Lhh_wealth_asinh            2539
## Lvill_eligible_ratio        2270
## Lsavingsstock_asinh3        3187
## Lconsumpti_x_Ldietarydi     2679
## Lconsumpti_x_Lproductiv     3445
## Ldietarydi_x_Lassetscon     3143
## 
## Family Specific Parameters: 
##       Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sigma     1.58      0.02     1.55     1.62 1.00     5780     3160
## 
## Draws were sampled using sampling(NUTS). For each parameter, Bulk_ESS
## and Tail_ESS are effective sample size measures, and Rhat is the potential
## scale reduction factor on split chains (at convergence, Rhat = 1).

Prior summery - how informative are priors

prior_summary(dietary_diversity_normal_bayesmodel)
##                 prior     class                    coef group resp dpar nlpar
##           normal(0,6)         b                                              
##           normal(0,6)         b          cost_deviation                      
##           normal(0,6)         b     dietarydiversity_R1                      
##           normal(0,6)         b Lconsumpti_x_Ldietarydi                      
##           normal(0,6)         b Lconsumpti_x_Lproductiv                      
##           normal(0,6)         b Ldietarydi_x_Lassetscon                      
##           normal(0,6)         b        Lhh_wealth_asinh                      
##           normal(0,6)         b    Lsavingsstock_asinh3                      
##           normal(0,6)         b    Lvill_eligible_ratio                      
##           normal(0,6)         b               treat_any                      
##           normal(0,6)         b                treat_GK                      
##  student_t(3, 5, 2.5) Intercept                                              
##  student_t(3, 0, 2.5)        sd                                              
##  student_t(3, 0, 2.5)        sd                         block                
##  student_t(3, 0, 2.5)        sd               Intercept block                
##  student_t(3, 0, 2.5)        sd                           vid                
##  student_t(3, 0, 2.5)        sd               Intercept   vid                
##  student_t(3, 0, 2.5)     sigma                                              
##  bound       source
##                user
##        (vectorized)
##        (vectorized)
##        (vectorized)
##        (vectorized)
##        (vectorized)
##        (vectorized)
##        (vectorized)
##        (vectorized)
##        (vectorized)
##        (vectorized)
##             default
##             default
##        (vectorized)
##        (vectorized)
##        (vectorized)
##        (vectorized)
##             default
check_prior(dietary_diversity_normal_bayesmodel)
##                    Parameter Prior_Quality
## 1                b_Intercept uninformative
## 2           b_cost_deviation uninformative
## 3                b_treat_any uninformative
## 4                 b_treat_GK uninformative
## 5      b_dietarydiversity_R1 uninformative
## 6         b_Lhh_wealth_asinh uninformative
## 7     b_Lvill_eligible_ratio uninformative
## 8     b_Lsavingsstock_asinh3 uninformative
## 9  b_Lconsumpti_x_Ldietarydi uninformative
## 10 b_Lconsumpti_x_Lproductiv uninformative
## 11 b_Ldietarydi_x_Lassetscon uninformative

Diagnostics

# trace diagnostic plot
mcmc_trace(dietary_diversity_normal_bayesmodel, n_warmup = 0,
           pars = c("b_Intercept", "b_cost_deviation", "b_treat_any", 
                    "b_treat_GK", "b_dietarydiversity_R1", "b_Lhh_wealth_asinh",
                    "b_Lvill_eligible_ratio", "b_Lsavingsstock_asinh3",
                    "b_Lconsumpti_x_Ldietarydi", "b_Lconsumpti_x_Lproductiv", 
                    "b_Ldietarydi_x_Lassetscon", "sd_block__Intercept", 
                    "sd_vid__Intercept", "sigma"))

ggsave("table_3_diagnostics\\dietary_div_normal_trace.png", plot = last_plot(), width = 12, height = 5)

#density diagnostic plot
mcmc_dens(dietary_diversity_normal_bayesmodel,
           pars = c("b_Intercept", "b_cost_deviation", "b_treat_any", 
                    "b_treat_GK", "b_dietarydiversity_R1", "b_Lhh_wealth_asinh",
                    "b_Lvill_eligible_ratio", "b_Lsavingsstock_asinh3",
                    "b_Lconsumpti_x_Ldietarydi", "b_Lconsumpti_x_Lproductiv", 
                    "b_Ldietarydi_x_Lassetscon", "sd_block__Intercept", 
                    "sd_vid__Intercept", "sigma"))

ggsave("table_3_diagnostics\\dietary_div_normal_dens.png", plot = last_plot(), width = 12, height = 5)

mcmc_dens_overlay(dietary_diversity_normal_bayesmodel,
           pars = c("b_Intercept", "b_cost_deviation", "b_treat_any", 
                    "b_treat_GK", "b_dietarydiversity_R1", "b_Lhh_wealth_asinh",
                    "b_Lvill_eligible_ratio", "b_Lsavingsstock_asinh3",
                    "b_Lconsumpti_x_Ldietarydi", "b_Lconsumpti_x_Lproductiv", 
                    "b_Ldietarydi_x_Lassetscon", "sd_block__Intercept", 
                    "sd_vid__Intercept", "sigma"))

ggsave("table_3_diagnostics\\dietary_div_normal_overlay.png", plot = last_plot(), width = 12, height = 5)


#acf (auto-correlation) diagnostic plot
mcmc_acf(dietary_diversity_normal_bayesmodel,
           pars = c("b_Intercept", "b_cost_deviation", "b_treat_any", 
                    "b_treat_GK", "b_dietarydiversity_R1", "b_Lhh_wealth_asinh",
                    "b_Lvill_eligible_ratio", "b_Lsavingsstock_asinh3",
                    "b_Lconsumpti_x_Ldietarydi", "b_Lconsumpti_x_Lproductiv", 
                    "b_Ldietarydi_x_Lassetscon", "sd_block__Intercept", 
                    "sd_vid__Intercept", "sigma"))

ggsave("table_3_diagnostics\\dietary_div_normal_acf.png", plot = last_plot(), width = 12, height = 5)

posterior predictive checks

pp_check(dietary_diversity_normal_bayesmodel, nsamples = 100)
## Warning: Argument 'nsamples' is deprecated. Please use argument 'ndraws'
## instead.

pp_check(dietary_diversity_normal_bayesmodel, nsamples = 10, type = 'error_scatter_avg', alpha = .1)
## Warning: Argument 'nsamples' is deprecated. Please use argument 'ndraws'
## instead.

Total Household Wealth - Primary Outcome

Load the total household wealth model variables

total household wealth: This is the basic benchmarking model utilizing brm() default, uninformed priors. Removed Lhh_wealth_asinh to account for collinearity issues.

hh_wealth_normal_bayesmodel <-
  brm(formula = wealth_asinh | weights(samp_wgt) ~
        cost_deviation + treat_any + treat_GK +
        wealth_asinh_R1 + Lvill_eligible_ratio + Lowndwelling +
        (1 | block) + (1 | vid),
     data = hh_wealth_data,
     family = gaussian("identity"),
     set_prior("normal(0,6)", class = "b"),
     seed = 1272022,
     warmup = 1000,
     iter = 2000,
     thin = 1,
     control = list(adapt_delta = .99, max_treedepth = 10),
     #backend = "cmdstanr",
     cores = 4, #overrides default 1 core
     #threads = 3,need to get cmdstanr package working here
     save_pars = save_pars(all = TRUE), # potentially allows for more post-processing functionality
     file = "informed_prior_outcomes\\hh_wealth_normal_bayes")

Model Summery

summary(hh_wealth_normal_bayesmodel)
##  Family: gaussian 
##   Links: mu = identity; sigma = identity 
## Formula: wealth_asinh | weights(samp_wgt) ~ cost_deviation + treat_any + treat_GK + wealth_asinh_R1 + Lvill_eligible_ratio + Lowndwelling + (1 | block) + (1 | vid) 
##    Data: hh_wealth_data (Number of observations: 1751) 
##   Draws: 4 chains, each with iter = 2000; warmup = 1000; thin = 1;
##          total post-warmup draws = 4000
## 
## Group-Level Effects: 
## ~block (Number of levels: 22) 
##               Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sd(Intercept)     0.52      0.15     0.24     0.84 1.00      838     1059
## 
## ~vid (Number of levels: 248) 
##               Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sd(Intercept)     1.03      0.09     0.85     1.22 1.00     1546     2483
## 
## Population-Level Effects: 
##                      Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS
## Intercept                7.81      0.32     7.18     8.42 1.00     3221
## cost_deviation           0.00      0.00    -0.00     0.00 1.00     3080
## treat_any                0.06      0.23    -0.40     0.52 1.00     2187
## treat_GK                 0.01      0.23    -0.43     0.47 1.00     2311
## wealth_asinh_R1          0.18      0.02     0.15     0.21 1.00     5499
## Lvill_eligible_ratio    -0.10      0.84    -1.79     1.49 1.00     1882
## Lowndwelling             3.46      0.20     3.07     3.86 1.00     5247
##                      Tail_ESS
## Intercept                3327
## cost_deviation           3530
## treat_any                2541
## treat_GK                 2728
## wealth_asinh_R1          3585
## Lvill_eligible_ratio     2751
## Lowndwelling             3514
## 
## Family Specific Parameters: 
##       Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sigma     3.43      0.04     3.35     3.51 1.00     6564     2972
## 
## Draws were sampled using sampling(NUTS). For each parameter, Bulk_ESS
## and Tail_ESS are effective sample size measures, and Rhat is the potential
## scale reduction factor on split chains (at convergence, Rhat = 1).

Prior summery - how informative are priors

prior_summary(hh_wealth_normal_bayesmodel)
##                    prior     class                 coef group resp dpar nlpar
##              normal(0,6)         b                                           
##              normal(0,6)         b       cost_deviation                      
##              normal(0,6)         b         Lowndwelling                      
##              normal(0,6)         b Lvill_eligible_ratio                      
##              normal(0,6)         b            treat_any                      
##              normal(0,6)         b             treat_GK                      
##              normal(0,6)         b      wealth_asinh_R1                      
##  student_t(3, 13.9, 2.5) Intercept                                           
##     student_t(3, 0, 2.5)        sd                                           
##     student_t(3, 0, 2.5)        sd                      block                
##     student_t(3, 0, 2.5)        sd            Intercept block                
##     student_t(3, 0, 2.5)        sd                        vid                
##     student_t(3, 0, 2.5)        sd            Intercept   vid                
##     student_t(3, 0, 2.5)     sigma                                           
##  bound       source
##                user
##        (vectorized)
##        (vectorized)
##        (vectorized)
##        (vectorized)
##        (vectorized)
##        (vectorized)
##             default
##             default
##        (vectorized)
##        (vectorized)
##        (vectorized)
##        (vectorized)
##             default
check_prior(hh_wealth_normal_bayesmodel)
##                Parameter Prior_Quality
## 1            b_Intercept uninformative
## 2       b_cost_deviation uninformative
## 3            b_treat_any uninformative
## 4             b_treat_GK uninformative
## 5      b_wealth_asinh_R1 uninformative
## 6 b_Lvill_eligible_ratio   informative
## 7         b_Lowndwelling uninformative

Diagnostics

# trace diagnostic plot
mcmc_trace(hh_wealth_normal_bayesmodel, n_warmup = 0,
           pars = c("b_Intercept", "b_cost_deviation", "b_treat_any", 
                    "b_treat_GK", "b_wealth_asinh_R1",
                    "b_Lvill_eligible_ratio", "b_Lowndwelling",
                    "sd_block__Intercept", "sd_vid__Intercept", "sigma"))

ggsave("table_3_diagnostics\\hh_wealth_normal_trace.png", plot = last_plot(), width = 12, height = 5)

#density diagnostic plot
mcmc_dens(hh_wealth_normal_bayesmodel,
           pars = c("b_Intercept", "b_cost_deviation", "b_treat_any", 
                    "b_treat_GK", "b_wealth_asinh_R1",
                    "b_Lvill_eligible_ratio", "b_Lowndwelling",
                    "sd_block__Intercept", "sd_vid__Intercept", "sigma"))

ggsave("table_3_diagnostics\\hh_wealth_normal_dens.png", plot = last_plot(), width = 12, height = 5)

mcmc_dens_overlay(hh_wealth_normal_bayesmodel,
           pars = c("b_Intercept", "b_cost_deviation", "b_treat_any", 
                    "b_treat_GK", "b_wealth_asinh_R1",
                    "b_Lvill_eligible_ratio", "b_Lowndwelling",
                    "sd_block__Intercept", "sd_vid__Intercept", "sigma"))

ggsave("table_3_diagnostics\\hh_wealth_normal_dens_overlay.png", plot = last_plot(), width = 12, height = 5)


#acf (auto-correlation) diagnostic plot
mcmc_acf(hh_wealth_normal_bayesmodel,
           pars = c("b_Intercept", "b_cost_deviation", "b_treat_any", 
                    "b_treat_GK", "b_wealth_asinh_R1",
                    "b_Lvill_eligible_ratio", "b_Lowndwelling",
                    "sd_block__Intercept", "sd_vid__Intercept", "sigma"))

ggsave("table_3_diagnostics\\hh_wealth_normal_acf.png", plot = last_plot(), width = 12, height = 5)

posterior predictive checks

pp_check(hh_wealth_normal_bayesmodel, nsamples = 100)
## Warning: Argument 'nsamples' is deprecated. Please use argument 'ndraws'
## instead.

pp_check(hh_wealth_normal_bayesmodel, nsamples = 10, type = 'error_scatter_avg', alpha = .1)
## Warning: Argument 'nsamples' is deprecated. Please use argument 'ndraws'
## instead.

Borrowing Stock - Secondary Outcome

Load the borrowing stock model variables

borrowing stock: This is the basic benchmarking model utilizing brm() global weakly informed priors with mean = 0 and 6 standard deviations.

borrowing_stock_normal_bayesmodel <-
  brm(formula = borrowingstock_asinh | weights(samp_wgt) ~
        cost_deviation + treat_any + treat_GK +
        borrowingstock_asinh_R1 + Lhh_wealth_asinh + Lvill_eligible_ratio +
        (1 | block) + (1 | vid),
     data = borrowing_stock_data,
     family = gaussian("identity"),
     set_prior("normal(0,6)", class = "b"),
     seed = 1272022,
     warmup = 1000,
     iter = 2000,
     thin = 1,
     control = list(adapt_delta = .95, max_treedepth = 10),
     #backend = "cmdstanr",
     cores = 4, #overrides default 1 core
     #threads = 3,need to get cmdstanr package working here
     save_pars = save_pars(all = TRUE), # potentially allows for more post-processing functionality
     file = "informed_prior_outcomes\\borrowing_stock_normal_bayes")

Model Summery

summary(borrowing_stock_normal_bayesmodel)
##  Family: gaussian 
##   Links: mu = identity; sigma = identity 
## Formula: borrowingstock_asinh | weights(samp_wgt) ~ cost_deviation + treat_any + treat_GK + borrowingstock_asinh_R1 + Lhh_wealth_asinh + Lvill_eligible_ratio + (1 | block) + (1 | vid) 
##    Data: borrowing_stock_data (Number of observations: 1751) 
##   Draws: 4 chains, each with iter = 2000; warmup = 1000; thin = 1;
##          total post-warmup draws = 4000
## 
## Group-Level Effects: 
## ~block (Number of levels: 22) 
##               Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sd(Intercept)     0.93      0.20     0.58     1.37 1.00     1115     1692
## 
## ~vid (Number of levels: 248) 
##               Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sd(Intercept)     1.19      0.11     0.99     1.41 1.00     1551     2673
## 
## Population-Level Effects: 
##                         Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS
## Intercept                   6.13      0.43     5.28     6.97 1.00     2114
## cost_deviation              0.00      0.00    -0.00     0.00 1.00     3212
## treat_any                  -0.73      0.28    -1.28    -0.19 1.00     2114
## treat_GK                    0.67      0.28     0.13     1.23 1.00     1993
## borrowingstock_asinh_R1     0.23      0.02     0.20     0.26 1.00     4909
## Lhh_wealth_asinh           -0.00      0.02    -0.04     0.03 1.00     4630
## Lvill_eligible_ratio       -0.50      1.08    -2.56     1.69 1.00     1775
##                         Tail_ESS
## Intercept                   2566
## cost_deviation              3023
## treat_any                   2525
## treat_GK                    2440
## borrowingstock_asinh_R1     3375
## Lhh_wealth_asinh            3116
## Lvill_eligible_ratio        2416
## 
## Family Specific Parameters: 
##       Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sigma     4.48      0.05     4.38     4.59 1.00     6601     2872
## 
## Draws were sampled using sampling(NUTS). For each parameter, Bulk_ESS
## and Tail_ESS are effective sample size measures, and Rhat is the potential
## scale reduction factor on split chains (at convergence, Rhat = 1).

Prior summery - how informative are priors

prior_summary(borrowing_stock_normal_bayesmodel)
##                   prior     class                    coef group resp dpar nlpar
##             normal(0,6)         b                                              
##             normal(0,6)         b borrowingstock_asinh_R1                      
##             normal(0,6)         b          cost_deviation                      
##             normal(0,6)         b        Lhh_wealth_asinh                      
##             normal(0,6)         b    Lvill_eligible_ratio                      
##             normal(0,6)         b               treat_any                      
##             normal(0,6)         b                treat_GK                      
##  student_t(3, 9.2, 2.9) Intercept                                              
##    student_t(3, 0, 2.9)        sd                                              
##    student_t(3, 0, 2.9)        sd                         block                
##    student_t(3, 0, 2.9)        sd               Intercept block                
##    student_t(3, 0, 2.9)        sd                           vid                
##    student_t(3, 0, 2.9)        sd               Intercept   vid                
##    student_t(3, 0, 2.9)     sigma                                              
##  bound       source
##                user
##        (vectorized)
##        (vectorized)
##        (vectorized)
##        (vectorized)
##        (vectorized)
##        (vectorized)
##             default
##             default
##        (vectorized)
##        (vectorized)
##        (vectorized)
##        (vectorized)
##             default
check_prior(borrowing_stock_normal_bayesmodel)
##                   Parameter Prior_Quality
## 1               b_Intercept uninformative
## 2          b_cost_deviation uninformative
## 3               b_treat_any uninformative
## 4                b_treat_GK uninformative
## 5 b_borrowingstock_asinh_R1 uninformative
## 6        b_Lhh_wealth_asinh uninformative
## 7    b_Lvill_eligible_ratio   informative

Diagnostics

# trace diagnostic plot
mcmc_trace(borrowing_stock_normal_bayesmodel, n_warmup = 0,
           pars = c("b_Intercept", "b_cost_deviation", "b_treat_any", "b_treat_GK",
                    "b_borrowingstock_asinh_R1", "b_Lhh_wealth_asinh", "b_Lvill_eligible_ratio",
                    "sd_block__Intercept", "sd_vid__Intercept", "sigma"))

ggsave("table_3_diagnostics\\borrowing_stock_normal_trace.png", plot = last_plot(), width = 12, height = 5)

#density diagnostic plots
mcmc_dens(borrowing_stock_normal_bayesmodel,
           pars = c("b_Intercept", "b_cost_deviation", "b_treat_any", "b_treat_GK",
                    "b_borrowingstock_asinh_R1", "b_Lhh_wealth_asinh", "b_Lvill_eligible_ratio",
                    "sd_block__Intercept", "sd_vid__Intercept", "sigma"))

ggsave("table_3_diagnostics\\borrowing_stock_normal_dens.png", plot = last_plot(), width = 12, height = 5)

mcmc_dens_overlay(borrowing_stock_normal_bayesmodel,
           pars = c("b_Intercept", "b_cost_deviation", "b_treat_any", "b_treat_GK",
                    "b_borrowingstock_asinh_R1", "b_Lhh_wealth_asinh", "b_Lvill_eligible_ratio",
                    "sd_block__Intercept", "sd_vid__Intercept", "sigma"))

ggsave("table_3_diagnostics\\borrowing_stock_normal_dens_overlay.png", plot = last_plot(), width = 12, height = 5)

#acf (auto-correlation) diagnostic plot
mcmc_acf(borrowing_stock_normal_bayesmodel,
           pars = c("b_Intercept", "b_cost_deviation", "b_treat_any", "b_treat_GK",
                    "b_borrowingstock_asinh_R1", "b_Lhh_wealth_asinh", "b_Lvill_eligible_ratio",
                    "sd_block__Intercept", "sd_vid__Intercept", "sigma"))

ggsave("table_3_diagnostics\\borrowing_stock_normal_acf.png", plot = last_plot(), width = 12, height = 5)

posterior predictive checks

pp_check(borrowing_stock_normal_bayesmodel, nsamples = 100)
## Warning: Argument 'nsamples' is deprecated. Please use argument 'ndraws'
## instead.

pp_check(borrowing_stock_normal_bayesmodel, nsamples = 10, type = 'error_scatter_avg', alpha = .1)
## Warning: Argument 'nsamples' is deprecated. Please use argument 'ndraws'
## instead.

Savings Stock - Secondary Outcome

Load the savings stock model variables

savings stock: This is the basic bechmarking model utilzing brm() default uninformed priors

savings_stock_normal_bayesmodel <-
  brm(formula = savingsstock_asinh | weights(samp_wgt) ~
        cost_deviation + treat_any + treat_GK +
        savingsstock_asinh_R1 + Lhh_wealth_asinh + Lvill_eligible_ratio + 
        Lconsumpti_x_Lproductiv + Lconsumpti_x_Lassetscon +
        (1 | block) + (1 | vid),
     data = savings_stock_data,
     family = gaussian("identity"),
     set_prior("normal(0,6)", class = "b"),
     seed = 1272022,
     warmup = 1000,
     iter = 2000,
     thin = 1,
     control = list(adapt_delta = .95, max_treedepth = 10),
     #backend = "cmdstanr",
     cores = 4, #overrides default 1 core
     #threads = 3,need to get cmdstanr package working here
     save_pars = save_pars(all = TRUE), # potentially allows for more post-processing functionality
     file = "informed_prior_outcomes\\savings_stock_normal_bayes")

Model Summery

summary(savings_stock_normal_bayesmodel)
##  Family: gaussian 
##   Links: mu = identity; sigma = identity 
## Formula: savingsstock_asinh | weights(samp_wgt) ~ cost_deviation + treat_any + treat_GK + savingsstock_asinh_R1 + Lhh_wealth_asinh + Lvill_eligible_ratio + Lconsumpti_x_Lproductiv + Lconsumpti_x_Lassetscon + (1 | block) + (1 | vid) 
##    Data: savings_stock_data (Number of observations: 1751) 
##   Draws: 4 chains, each with iter = 2000; warmup = 1000; thin = 1;
##          total post-warmup draws = 4000
## 
## Group-Level Effects: 
## ~block (Number of levels: 22) 
##               Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sd(Intercept)     0.70      0.19     0.36     1.12 1.00      840     1205
## 
## ~vid (Number of levels: 248) 
##               Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sd(Intercept)     1.40      0.11     1.18     1.63 1.01     1074     1983
## 
## Population-Level Effects: 
##                         Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS
## Intercept                   1.78      0.46     0.88     2.68 1.00     2347
## cost_deviation              0.00      0.00    -0.00     0.00 1.00     2022
## treat_any                  -0.13      0.30    -0.73     0.48 1.00     1409
## treat_GK                    1.31      0.30     0.73     1.92 1.00     1287
## savingsstock_asinh_R1       0.17      0.01     0.14     0.20 1.00     4002
## Lhh_wealth_asinh           -0.05      0.02    -0.09    -0.02 1.00     4623
## Lvill_eligible_ratio        1.54      1.17    -0.68     3.89 1.00     1232
## Lconsumpti_x_Lproductiv     0.02      0.00     0.02     0.03 1.00     5177
## Lconsumpti_x_Lassetscon     0.01      0.00     0.01     0.01 1.00     5780
##                         Tail_ESS
## Intercept                   3166
## cost_deviation              2700
## treat_any                   1854
## treat_GK                    2084
## savingsstock_asinh_R1       2628
## Lhh_wealth_asinh            3062
## Lvill_eligible_ratio        1905
## Lconsumpti_x_Lproductiv     3323
## Lconsumpti_x_Lassetscon     3550
## 
## Family Specific Parameters: 
##       Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sigma     4.15      0.05     4.06     4.25 1.00     5102     3097
## 
## Draws were sampled using sampling(NUTS). For each parameter, Bulk_ESS
## and Tail_ESS are effective sample size measures, and Rhat is the potential
## scale reduction factor on split chains (at convergence, Rhat = 1).

Prior summery - how informative are priors

prior_summary(savings_stock_normal_bayesmodel)
##                   prior     class                    coef group resp dpar nlpar
##             normal(0,6)         b                                              
##             normal(0,6)         b          cost_deviation                      
##             normal(0,6)         b Lconsumpti_x_Lassetscon                      
##             normal(0,6)         b Lconsumpti_x_Lproductiv                      
##             normal(0,6)         b        Lhh_wealth_asinh                      
##             normal(0,6)         b    Lvill_eligible_ratio                      
##             normal(0,6)         b   savingsstock_asinh_R1                      
##             normal(0,6)         b               treat_any                      
##             normal(0,6)         b                treat_GK                      
##  student_t(3, 8.7, 3.1) Intercept                                              
##    student_t(3, 0, 3.1)        sd                                              
##    student_t(3, 0, 3.1)        sd                         block                
##    student_t(3, 0, 3.1)        sd               Intercept block                
##    student_t(3, 0, 3.1)        sd                           vid                
##    student_t(3, 0, 3.1)        sd               Intercept   vid                
##    student_t(3, 0, 3.1)     sigma                                              
##  bound       source
##                user
##        (vectorized)
##        (vectorized)
##        (vectorized)
##        (vectorized)
##        (vectorized)
##        (vectorized)
##        (vectorized)
##        (vectorized)
##             default
##             default
##        (vectorized)
##        (vectorized)
##        (vectorized)
##        (vectorized)
##             default
check_prior(savings_stock_normal_bayesmodel)
##                   Parameter Prior_Quality
## 1               b_Intercept uninformative
## 2          b_cost_deviation uninformative
## 3               b_treat_any uninformative
## 4                b_treat_GK uninformative
## 5   b_savingsstock_asinh_R1 uninformative
## 6        b_Lhh_wealth_asinh uninformative
## 7    b_Lvill_eligible_ratio   informative
## 8 b_Lconsumpti_x_Lproductiv uninformative
## 9 b_Lconsumpti_x_Lassetscon uninformative

Diagnostics

# trace diagnostic plot
mcmc_trace(savings_stock_normal_bayesmodel, n_warmup = 0,
           pars = c("b_Intercept", "b_cost_deviation", "b_treat_any", "b_treat_GK",
                    "b_savingsstock_asinh_R1", "b_Lhh_wealth_asinh", "b_Lvill_eligible_ratio",
                    "b_Lconsumpti_x_Lproductiv", "b_Lconsumpti_x_Lassetscon",
                    "sd_block__Intercept", "sd_vid__Intercept", "sigma"))

ggsave("table_3_diagnostics\\savings_stock_normal_trace.png", plot = last_plot(), width = 12, height = 5)

#density diagnostic plots
mcmc_dens(savings_stock_normal_bayesmodel,
           pars = c("b_Intercept", "b_cost_deviation", "b_treat_any", "b_treat_GK",
                    "b_savingsstock_asinh_R1", "b_Lhh_wealth_asinh", "b_Lvill_eligible_ratio",
                    "b_Lconsumpti_x_Lproductiv", "b_Lconsumpti_x_Lassetscon",
                    "sd_block__Intercept", "sd_vid__Intercept", "sigma"))

ggsave("table_3_diagnostics\\savings_stock_normal_dens.png", plot = last_plot(), width = 12, height = 5)

mcmc_dens_overlay(savings_stock_normal_bayesmodel,
           pars = c("b_Intercept", "b_cost_deviation", "b_treat_any", "b_treat_GK",
                    "b_savingsstock_asinh_R1", "b_Lhh_wealth_asinh", "b_Lvill_eligible_ratio",
                    "b_Lconsumpti_x_Lproductiv", "b_Lconsumpti_x_Lassetscon",
                    "sd_block__Intercept", "sd_vid__Intercept", "sigma"))

ggsave("table_3_diagnostics\\savings_stock_normal_dens_overlay.png", plot = last_plot(), width = 12, height = 5)

#acf (auto-correlation) diagnostic plot
mcmc_acf(savings_stock_normal_bayesmodel,
           pars = c("b_Intercept", "b_cost_deviation", "b_treat_any", "b_treat_GK",
                    "b_savingsstock_asinh_R1", "b_Lhh_wealth_asinh", "b_Lvill_eligible_ratio",
                    "b_Lconsumpti_x_Lproductiv", "b_Lconsumpti_x_Lassetscon",
                    "sd_block__Intercept", "sd_vid__Intercept", "sigma"))

ggsave("table_3_diagnostics\\savings_stock_normal_acf.png", plot = last_plot(), width = 12, height = 5)

posterior predictive checks

pp_check(savings_stock_normal_bayesmodel, nsamples = 100)
## Warning: Argument 'nsamples' is deprecated. Please use argument 'ndraws'
## instead.

pp_check(savings_stock_normal_bayesmodel, nsamples = 10, type = 'error_scatter_avg', alpha = .1)
## Warning: Argument 'nsamples' is deprecated. Please use argument 'ndraws'
## instead.

Health Knowledge - Secondary Outcome

Load the health knowledge model variables

health knowledge: This is the basic bechmarking model utilzing brm() default, uninformed priors

health_knowledge_normal_bayesmodel <-
  brm(formula = health_knowledge | weights(samp_wgt) ~
        cost_deviation + treat_any + treat_GK +
        health_knowledge_R1 + 
        Lhh_wealth_asinh + Lvill_eligible_ratio +
        (1 | block) + (1 | vid),
     data = health_knowledge_data,
     family = gaussian("identity"),
     set_prior("normal(0,6)", class = "b"),
     seed = 1272022,
     warmup = 1000,
     iter = 2000,
     thin = 1,
     control = list(adapt_delta = .95, max_treedepth = 10),
     #backend = "cmdstanr",
     cores = 4, #overrides default 1 core
     #threads = 3,need to get cmdstanr package working here
     save_pars = save_pars(all = TRUE), # potentially allows for more post-processing functionality
     file = "informed_prior_outcomes\\health_knowledge_normal_bayes")

Model Summery

summary(health_knowledge_normal_bayesmodel)
##  Family: gaussian 
##   Links: mu = identity; sigma = identity 
## Formula: health_knowledge | weights(samp_wgt) ~ cost_deviation + treat_any + treat_GK + health_knowledge_R1 + Lhh_wealth_asinh + Lvill_eligible_ratio + (1 | block) + (1 | vid) 
##    Data: health_knowledge_data (Number of observations: 1751) 
##   Draws: 4 chains, each with iter = 2000; warmup = 1000; thin = 1;
##          total post-warmup draws = 4000
## 
## Group-Level Effects: 
## ~block (Number of levels: 22) 
##               Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sd(Intercept)     0.34      0.22     0.02     0.81 1.00      454     1098
## 
## ~vid (Number of levels: 248) 
##               Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sd(Intercept)     1.64      0.12     1.42     1.88 1.00     1474     2879
## 
## Population-Level Effects: 
##                      Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS
## Intercept                1.37      0.40     0.60     2.17 1.00     2668
## cost_deviation           0.00      0.00    -0.00     0.00 1.00     2071
## treat_any                0.10      0.33    -0.54     0.76 1.00     1447
## treat_GK                -0.08      0.34    -0.76     0.60 1.00     1364
## health_knowledge_R1      0.06      0.02     0.03     0.09 1.00     4833
## Lhh_wealth_asinh         0.11      0.02     0.07     0.15 1.00     5603
## Lvill_eligible_ratio     0.61      1.15    -1.66     2.89 1.00     1349
##                      Tail_ESS
## Intercept                2726
## cost_deviation           2465
## treat_any                2351
## treat_GK                 2339
## health_knowledge_R1      2690
## Lhh_wealth_asinh         2820
## Lvill_eligible_ratio     1861
## 
## Family Specific Parameters: 
##       Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sigma     4.35      0.05     4.25     4.45 1.00     4968     2742
## 
## Draws were sampled using sampling(NUTS). For each parameter, Bulk_ESS
## and Tail_ESS are effective sample size measures, and Rhat is the potential
## scale reduction factor on split chains (at convergence, Rhat = 1).

Prior summery - how informative are priors

prior_summary(health_knowledge_normal_bayesmodel)
##                   prior     class                 coef group resp dpar nlpar
##             normal(0,6)         b                                           
##             normal(0,6)         b       cost_deviation                      
##             normal(0,6)         b  health_knowledge_R1                      
##             normal(0,6)         b     Lhh_wealth_asinh                      
##             normal(0,6)         b Lvill_eligible_ratio                      
##             normal(0,6)         b            treat_any                      
##             normal(0,6)         b             treat_GK                      
##  student_t(3, 3.2, 3.8) Intercept                                           
##    student_t(3, 0, 3.8)        sd                                           
##    student_t(3, 0, 3.8)        sd                      block                
##    student_t(3, 0, 3.8)        sd            Intercept block                
##    student_t(3, 0, 3.8)        sd                        vid                
##    student_t(3, 0, 3.8)        sd            Intercept   vid                
##    student_t(3, 0, 3.8)     sigma                                           
##  bound       source
##                user
##        (vectorized)
##        (vectorized)
##        (vectorized)
##        (vectorized)
##        (vectorized)
##        (vectorized)
##             default
##             default
##        (vectorized)
##        (vectorized)
##        (vectorized)
##        (vectorized)
##             default
check_prior(health_knowledge_normal_bayesmodel)
##                Parameter Prior_Quality
## 1            b_Intercept   informative
## 2       b_cost_deviation uninformative
## 3            b_treat_any uninformative
## 4             b_treat_GK uninformative
## 5  b_health_knowledge_R1 uninformative
## 6     b_Lhh_wealth_asinh uninformative
## 7 b_Lvill_eligible_ratio   informative

Diagnostics

# trace diagnostic plot
mcmc_trace(health_knowledge_normal_bayesmodel, n_warmup = 0,
           pars = c("b_Intercept", "b_cost_deviation", "b_treat_any", "b_treat_GK",
                    "b_health_knowledge_R1", "b_Lhh_wealth_asinh", "b_Lvill_eligible_ratio",
                    "sd_block__Intercept", "sd_vid__Intercept", "sigma"))

ggsave("table_3_diagnostics\\health_knowledge_normal_trace.png", plot = last_plot(), width = 12, height = 5)

#density diagnostic plots
mcmc_dens(health_knowledge_normal_bayesmodel,
           pars = c("b_Intercept", "b_cost_deviation", "b_treat_any", "b_treat_GK",
                    "b_health_knowledge_R1", "b_Lhh_wealth_asinh", "b_Lvill_eligible_ratio",
                    "sd_block__Intercept", "sd_vid__Intercept", "sigma"))

ggsave("table_3_diagnostics\\health_knowledge_normal_dens.png", plot = last_plot(), width = 12, height = 5)

mcmc_dens_overlay(health_knowledge_normal_bayesmodel,
           pars = c("b_Intercept", "b_cost_deviation", "b_treat_any", "b_treat_GK",
                    "b_health_knowledge_R1", "b_Lhh_wealth_asinh", "b_Lvill_eligible_ratio",
                    "sd_block__Intercept", "sd_vid__Intercept", "sigma"))

ggsave("table_3_diagnostics\\health_knowledge_normal_dens_overlay.png", plot = last_plot(), width = 12, height = 5)

#acf (auto-correlation) diagnostic plot
mcmc_acf(health_knowledge_normal_bayesmodel,
           pars = c("b_Intercept", "b_cost_deviation", "b_treat_any", "b_treat_GK",
                    "b_health_knowledge_R1", "b_Lhh_wealth_asinh", "b_Lvill_eligible_ratio",
                    "sd_block__Intercept", "sd_vid__Intercept", "sigma"))

ggsave("table_3_diagnostics\\health_knowledge_normal_acf.png", plot = last_plot(), width = 12, height = 5)

posterior predictive checks

pp_check(health_knowledge_normal_bayesmodel, ndraws = 100)

pp_check(health_knowledge_normal_bayesmodel, ndraws = 10, type = 'error_scatter_avg', alpha = .1)

sanitation practices - Secondary Outcome

Load the sanitation practices model variables

sanitation practices: This is the basic benchmarking model utilizing the default, uninformed priors

sanitation_practices_normal_bayesmodel <-
  brm(formula = sanitation_practices | weights(samp_wgt) ~
        cost_deviation + treat_any + treat_GK +
        sanitation_practices_R1 + Lhh_wealth_asinh + Lvill_eligible_ratio +
        Lproductiv_x_Lassetscon +
        (1 | block) + (1 | vid),
     data = sanitation_practices_data,
     family = gaussian("identity"),
     set_prior("normal(0,6)", class = "b"),
     seed = 1272022,
     warmup = 1000,
     iter = 2000,
     thin = 1,
     control = list(adapt_delta = .95, max_treedepth = 10),
     #backend = "cmdstanr",
     cores = 4, #overrides default 1 core
     #threads = 3,need to get cmdstanr package working here
     save_pars = save_pars(all = TRUE), # potentially allows for more post-processing functionality
     file = "informed_prior_outcomes\\sanitation_practices_normal_bayes")

Model Summery

summary(sanitation_practices_normal_bayesmodel)
##  Family: gaussian 
##   Links: mu = identity; sigma = identity 
## Formula: sanitation_practices | weights(samp_wgt) ~ cost_deviation + treat_any + treat_GK + sanitation_practices_R1 + Lhh_wealth_asinh + Lvill_eligible_ratio + Lproductiv_x_Lassetscon + (1 | block) + (1 | vid) 
##    Data: sanitation_practices_data (Number of observations: 1751) 
##   Draws: 4 chains, each with iter = 2000; warmup = 1000; thin = 1;
##          total post-warmup draws = 4000
## 
## Group-Level Effects: 
## ~block (Number of levels: 22) 
##               Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sd(Intercept)     0.18      0.12     0.01     0.44 1.00      695     1527
## 
## ~vid (Number of levels: 248) 
##               Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sd(Intercept)     0.97      0.07     0.84     1.10 1.00     1975     3150
## 
## Population-Level Effects: 
##                         Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS
## Intercept                  -2.01      0.24    -2.49    -1.54 1.00     2838
## cost_deviation              0.00      0.00    -0.00     0.00 1.00     2594
## treat_any                   0.21      0.20    -0.18     0.61 1.00     1824
## treat_GK                   -0.31      0.20    -0.70     0.09 1.00     1809
## sanitation_practices_R1     0.11      0.02     0.07     0.14 1.00     6618
## Lhh_wealth_asinh            0.04      0.01     0.01     0.06 1.00     8011
## Lvill_eligible_ratio        0.03      0.69    -1.32     1.41 1.00     1698
## Lproductiv_x_Lassetscon     0.01      0.00     0.01     0.01 1.00     5045
##                         Tail_ESS
## Intercept                   3352
## cost_deviation              2614
## treat_any                   2550
## treat_GK                    2471
## sanitation_practices_R1     3076
## Lhh_wealth_asinh            3119
## Lvill_eligible_ratio        2384
## Lproductiv_x_Lassetscon     3458
## 
## Family Specific Parameters: 
##       Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sigma     2.63      0.03     2.57     2.69 1.00     7096     2932
## 
## Draws were sampled using sampling(NUTS). For each parameter, Bulk_ESS
## and Tail_ESS are effective sample size measures, and Rhat is the potential
## scale reduction factor on split chains (at convergence, Rhat = 1).

Prior summery - how informative are priors

prior_summary(sanitation_practices_normal_bayesmodel)
##                 prior     class                    coef group resp dpar nlpar
##           normal(0,6)         b                                              
##           normal(0,6)         b          cost_deviation                      
##           normal(0,6)         b        Lhh_wealth_asinh                      
##           normal(0,6)         b Lproductiv_x_Lassetscon                      
##           normal(0,6)         b    Lvill_eligible_ratio                      
##           normal(0,6)         b sanitation_practices_R1                      
##           normal(0,6)         b               treat_any                      
##           normal(0,6)         b                treat_GK                      
##  student_t(3, 0.1, 3) Intercept                                              
##    student_t(3, 0, 3)        sd                                              
##    student_t(3, 0, 3)        sd                         block                
##    student_t(3, 0, 3)        sd               Intercept block                
##    student_t(3, 0, 3)        sd                           vid                
##    student_t(3, 0, 3)        sd               Intercept   vid                
##    student_t(3, 0, 3)     sigma                                              
##  bound       source
##                user
##        (vectorized)
##        (vectorized)
##        (vectorized)
##        (vectorized)
##        (vectorized)
##        (vectorized)
##        (vectorized)
##             default
##             default
##        (vectorized)
##        (vectorized)
##        (vectorized)
##        (vectorized)
##             default
check_prior(sanitation_practices_normal_bayesmodel)
##                   Parameter Prior_Quality
## 1               b_Intercept   informative
## 2          b_cost_deviation uninformative
## 3               b_treat_any uninformative
## 4                b_treat_GK uninformative
## 5 b_sanitation_practices_R1 uninformative
## 6        b_Lhh_wealth_asinh uninformative
## 7    b_Lvill_eligible_ratio   informative
## 8 b_Lproductiv_x_Lassetscon uninformative

Diagnostics

# trace diagnostic plot
mcmc_trace(sanitation_practices_normal_bayesmodel, n_warmup = 0,
           pars = c("b_Intercept", "b_cost_deviation", "b_treat_any", "b_treat_GK",
                    "b_sanitation_practices_R1", "b_Lhh_wealth_asinh", "b_Lvill_eligible_ratio",
                    "b_Lproductiv_x_Lassetscon",
                    "sd_block__Intercept", "sd_vid__Intercept", "sigma"))

ggsave("table_3_diagnostics\\sanitation_practices_normal_trace.png", plot = last_plot(), width = 12, height = 5)

#density diagnostic plots
mcmc_dens(sanitation_practices_normal_bayesmodel,
           pars = c("b_Intercept", "b_cost_deviation", "b_treat_any", "b_treat_GK",
                    "b_sanitation_practices_R1", "b_Lhh_wealth_asinh", "b_Lvill_eligible_ratio",
                    "b_Lproductiv_x_Lassetscon",
                    "sd_block__Intercept", "sd_vid__Intercept", "sigma"))

ggsave("table_3_diagnostics\\sanitation_practices_normal_dens.png", plot = last_plot(), width = 12, height = 5)

mcmc_dens_overlay(sanitation_practices_normal_bayesmodel,
           pars = c("b_Intercept", "b_cost_deviation", "b_treat_any", "b_treat_GK",
                    "b_sanitation_practices_R1", "b_Lhh_wealth_asinh", "b_Lvill_eligible_ratio",
                    "b_Lproductiv_x_Lassetscon",
                    "sd_block__Intercept", "sd_vid__Intercept", "sigma"))

ggsave("table_3_diagnostics\\sanitation_practices_normal_dens_overlay.png", plot = last_plot(), width = 12, height = 5)

#acf (auto-correlation) diagnostic plot
mcmc_acf(sanitation_practices_normal_bayesmodel,
           pars = c("b_Intercept", "b_cost_deviation", "b_treat_any", "b_treat_GK",
                    "b_sanitation_practices_R1", "b_Lhh_wealth_asinh", "b_Lvill_eligible_ratio",
                    "b_Lproductiv_x_Lassetscon",
                    "sd_block__Intercept", "sd_vid__Intercept", "sigma"))

ggsave("table_3_diagnostics\\sanitation_practices_normal_acf.png", plot = last_plot(), width = 12, height = 5)

posterior predictive checks

pp_check(sanitation_practices_normal_bayesmodel, ndraws = 100)

pp_check(sanitation_practices_normal_bayesmodel, ndraws = 10, type = 'error_scatter_avg', alpha = .1)

productive assets - Secondary Outcome

Load the productive assets model variables

productive assets: This is the basic benchmarking model utilizing the default, uninformed priors

productive_assets_normal_bayesmodel <-
  brm(formula = productiveassets_asinh | weights(samp_wgt) ~
        cost_deviation + treat_any + treat_GK +
        productiveassets_asinh_R1 + Lhh_wealth_asinh + Lvill_eligible_ratio +
        Lconsumpti_x_Lassetscon +
        (1 | block) + (1 | vid),
     data = productive_assets_data,
     family = gaussian("identity"),
     set_prior("normal(0,6)", class = "b"),
     seed = 1272022,
     warmup = 1000,
     iter = 2000,
     thin = 1,
     control = list(adapt_delta = .95, max_treedepth = 10),
     #backend = "cmdstanr",
     cores = 4, #overrides default 1 core
     #threads = 3,need to get cmdstanr package working here
     save_pars = save_pars(all = TRUE), # potentially allows for more post-processing functionality
     file = "informed_prior_outcomes\\productive_assets_normal_bayes")

Model Summery

summary(productive_assets_normal_bayesmodel)
##  Family: gaussian 
##   Links: mu = identity; sigma = identity 
## Formula: productiveassets_asinh | weights(samp_wgt) ~ cost_deviation + treat_any + treat_GK + productiveassets_asinh_R1 + Lhh_wealth_asinh + Lvill_eligible_ratio + Lconsumpti_x_Lassetscon + (1 | block) + (1 | vid) 
##    Data: productive_assets_data (Number of observations: 1751) 
##   Draws: 4 chains, each with iter = 2000; warmup = 1000; thin = 1;
##          total post-warmup draws = 4000
## 
## Group-Level Effects: 
## ~block (Number of levels: 22) 
##               Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sd(Intercept)     0.34      0.07     0.22     0.50 1.00     1918     2377
## 
## ~vid (Number of levels: 248) 
##               Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sd(Intercept)     0.43      0.04     0.36     0.51 1.00     1525     2705
## 
## Population-Level Effects: 
##                           Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS
## Intercept                     6.89      0.21     6.48     7.31 1.00     4710
## cost_deviation                0.00      0.00     0.00     0.00 1.00     4761
## treat_any                     0.37      0.10     0.18     0.57 1.00     3640
## treat_GK                     -0.27      0.10    -0.46    -0.07 1.00     3640
## productiveassets_asinh_R1     0.33      0.02     0.30     0.37 1.00     9267
## Lhh_wealth_asinh             -0.00      0.01    -0.01     0.01 1.00     8694
## Lvill_eligible_ratio          0.12      0.40    -0.67     0.88 1.00     3211
## Lconsumpti_x_Lassetscon       0.01      0.00     0.00     0.01 1.00     4417
##                           Tail_ESS
## Intercept                     3456
## cost_deviation                3691
## treat_any                     3388
## treat_GK                      3483
## productiveassets_asinh_R1     3190
## Lhh_wealth_asinh              3238
## Lvill_eligible_ratio          3084
## Lconsumpti_x_Lassetscon       3279
## 
## Family Specific Parameters: 
##       Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sigma     1.50      0.02     1.47     1.53 1.00     7656     3088
## 
## Draws were sampled using sampling(NUTS). For each parameter, Bulk_ESS
## and Tail_ESS are effective sample size measures, and Rhat is the potential
## scale reduction factor on split chains (at convergence, Rhat = 1).

Prior summery - how informative are priors

prior_summary(productive_assets_normal_bayesmodel)
##                    prior     class                      coef group resp dpar
##              normal(0,6)         b                                          
##              normal(0,6)         b            cost_deviation                
##              normal(0,6)         b   Lconsumpti_x_Lassetscon                
##              normal(0,6)         b          Lhh_wealth_asinh                
##              normal(0,6)         b      Lvill_eligible_ratio                
##              normal(0,6)         b productiveassets_asinh_R1                
##              normal(0,6)         b                 treat_any                
##              normal(0,6)         b                  treat_GK                
##  student_t(3, 11.4, 2.5) Intercept                                          
##     student_t(3, 0, 2.5)        sd                                          
##     student_t(3, 0, 2.5)        sd                           block          
##     student_t(3, 0, 2.5)        sd                 Intercept block          
##     student_t(3, 0, 2.5)        sd                             vid          
##     student_t(3, 0, 2.5)        sd                 Intercept   vid          
##     student_t(3, 0, 2.5)     sigma                                          
##  nlpar bound       source
##                      user
##              (vectorized)
##              (vectorized)
##              (vectorized)
##              (vectorized)
##              (vectorized)
##              (vectorized)
##              (vectorized)
##                   default
##                   default
##              (vectorized)
##              (vectorized)
##              (vectorized)
##              (vectorized)
##                   default
check_prior(productive_assets_normal_bayesmodel)
##                     Parameter Prior_Quality
## 1                 b_Intercept uninformative
## 2            b_cost_deviation uninformative
## 3                 b_treat_any uninformative
## 4                  b_treat_GK uninformative
## 5 b_productiveassets_asinh_R1 uninformative
## 6          b_Lhh_wealth_asinh uninformative
## 7      b_Lvill_eligible_ratio uninformative
## 8   b_Lconsumpti_x_Lassetscon uninformative

Diagnostics

# trace diagnostic plot
mcmc_trace(productive_assets_normal_bayesmodel, n_warmup = 0,
           pars = c("b_Intercept", "b_cost_deviation", "b_treat_any", "b_treat_GK",
                    "b_productiveassets_asinh_R1", "b_Lhh_wealth_asinh", "b_Lvill_eligible_ratio",
                    "b_Lconsumpti_x_Lassetscon",
                    "sd_block__Intercept", "sd_vid__Intercept", "sigma"))

ggsave("table_3_diagnostics\\productive_assets_normal_trace.png", plot = last_plot(), width = 12, height = 5)

#density diagnostic plots
mcmc_dens(productive_assets_normal_bayesmodel,
           pars = c("b_Intercept", "b_cost_deviation", "b_treat_any", "b_treat_GK",
                    "b_productiveassets_asinh_R1", "b_Lhh_wealth_asinh", "b_Lvill_eligible_ratio",
                    "b_Lconsumpti_x_Lassetscon",
                    "sd_block__Intercept", "sd_vid__Intercept", "sigma"))

ggsave("table_3_diagnostics\\productive_assets_normal_dens.png", plot = last_plot(), width = 12, height = 5)

mcmc_dens_overlay(productive_assets_normal_bayesmodel,
           pars = c("b_Intercept", "b_cost_deviation", "b_treat_any", "b_treat_GK",
                    "b_productiveassets_asinh_R1", "b_Lhh_wealth_asinh", "b_Lvill_eligible_ratio",
                    "b_Lconsumpti_x_Lassetscon",
                    "sd_block__Intercept", "sd_vid__Intercept", "sigma"))

ggsave("table_3_diagnostics\\productive_assets_normal_dens_overlay.png", plot = last_plot(), width = 12, height = 5)

#acf (auto-correlation) diagnostic plot
mcmc_acf(productive_assets_normal_bayesmodel,
           pars = c("b_Intercept", "b_cost_deviation", "b_treat_any", "b_treat_GK",
                    "b_productiveassets_asinh_R1", "b_Lhh_wealth_asinh", "b_Lvill_eligible_ratio",
                    "b_Lconsumpti_x_Lassetscon",
                    "sd_block__Intercept", "sd_vid__Intercept", "sigma"))

ggsave("table_3_diagnostics\\productive_assets_normal_cf.png", plot = last_plot(), width = 12, height = 5)

posterior predictive checks

pp_check(productive_assets_normal_bayesmodel, ndraws = 100)

pp_check(productive_assets_normal_bayesmodel, ndraws = 10, type = 'error_scatter_avg', alpha = .1)

Consumption Assets - Secondary Outcome

Load the consumption assets model variables

consumption assets: This is the basic benchmarking model utilizing a global normal prior with mean = and 6 standard deviations. This model has a number of errors: Warning:

consumption_assets_normal_bayesmodel <-
  brm(formula = assetsconsumption_asinh | weights(samp_wgt) ~
        cost_deviation + treat_any + treat_GK +
        assetsconsumption_asinh_R1 + Lhh_wealth_asinh + Lvill_eligible_ratio + Lroomsnumb + Ldurablesexpenditure +
        Ldietarydi_x_Lassetscon + Lproductiv_x_Lassetscon +
        (1 | block) + (1 | vid),
     data = consumption_assets_data,
     family = gaussian("identity"),
     set_prior("normal(0,6)", class = "b"),
     seed = 1272022,
     warmup = 1000,
     iter = 2000,
     thin = 1,
     control = list(adapt_delta = .95, max_treedepth = 10),
     #backend = "cmdstanr",
     cores = 4, #overrides default 1 core
     #threads = 3,need to get cmdstanr package working here
     save_pars = save_pars(all = TRUE), # potentially allows for more post-processing functionality
     file = "informed_prior_outcomes\\consumption_assets_normal_bayes")

Model Summery

summary(consumption_assets_normal_bayesmodel)
## Warning: Parts of the model have not converged (some Rhats are > 1.05). Be
## careful when analysing the results! We recommend running more iterations and/or
## setting stronger priors.
##  Family: gaussian 
##   Links: mu = identity; sigma = identity 
## Formula: assetsconsumption_asinh | weights(samp_wgt) ~ cost_deviation + treat_any + treat_GK + assetsconsumption_asinh_R1 + Lhh_wealth_asinh + Lvill_eligible_ratio + Lroomsnumb + Ldurablesexpenditure + Ldietarydi_x_Lassetscon + Lproductiv_x_Lassetscon + (1 | block) + (1 | vid) 
##    Data: consumption_assets_data (Number of observations: 1751) 
##   Draws: 4 chains, each with iter = 2000; warmup = 1000; thin = 1;
##          total post-warmup draws = 4000
## 
## Group-Level Effects: 
## ~block (Number of levels: 22) 
##               Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sd(Intercept)     0.89      1.09     0.00     2.72 3.49        4       14
## 
## ~vid (Number of levels: 248) 
##               Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sd(Intercept)     0.82      0.55     0.33     1.71 4.43        4       11
## 
## Population-Level Effects: 
##                            Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS
## Intercept                     -0.35      5.35   -15.70     4.06 3.42        4
## cost_deviation                 0.00      0.00     0.00     0.01 2.12        5
## treat_any                      0.54      1.22    -1.47     1.88 3.51        4
## treat_GK                       0.93      1.07    -1.23     1.79 3.56        4
## assetsconsumption_asinh_R1     0.34      0.35    -0.63     0.82 2.46        5
## Lhh_wealth_asinh              -0.17      0.23    -0.57    -0.01 1.60        7
## Lvill_eligible_ratio           0.70      0.50    -0.15     1.36 3.28        5
## Lroomsnumb                    -0.16      0.81    -1.49     0.95 2.93        5
## Ldurablesexpenditure          -0.00      0.00    -0.00     0.00 1.72        6
## Ldietarydi_x_Lassetscon        0.18      0.29     0.01     0.70 1.74        6
## Lproductiv_x_Lassetscon       -0.01      0.07    -0.18     0.08 2.72        5
##                            Tail_ESS
## Intercept                        11
## cost_deviation                   33
## treat_any                        11
## treat_GK                         13
## assetsconsumption_asinh_R1       11
## Lhh_wealth_asinh                 15
## Lvill_eligible_ratio             15
## Lroomsnumb                       11
## Ldurablesexpenditure             11
## Ldietarydi_x_Lassetscon          11
## Lproductiv_x_Lassetscon          11
## 
## Family Specific Parameters: 
##       Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sigma     7.29      6.27     3.33    19.92 2.53        6       13
## 
## Draws were sampled using sampling(NUTS). For each parameter, Bulk_ESS
## and Tail_ESS are effective sample size measures, and Rhat is the potential
## scale reduction factor on split chains (at convergence, Rhat = 1).

Prior summery - how informative are priors

prior_summary(consumption_assets_normal_bayesmodel)
##                   prior     class                       coef group resp dpar
##             normal(0,6)         b                                           
##             normal(0,6)         b assetsconsumption_asinh_R1                
##             normal(0,6)         b             cost_deviation                
##             normal(0,6)         b    Ldietarydi_x_Lassetscon                
##             normal(0,6)         b       Ldurablesexpenditure                
##             normal(0,6)         b           Lhh_wealth_asinh                
##             normal(0,6)         b    Lproductiv_x_Lassetscon                
##             normal(0,6)         b                 Lroomsnumb                
##             normal(0,6)         b       Lvill_eligible_ratio                
##             normal(0,6)         b                  treat_any                
##             normal(0,6)         b                   treat_GK                
##  student_t(3, 9.8, 2.5) Intercept                                           
##    student_t(3, 0, 2.5)        sd                                           
##    student_t(3, 0, 2.5)        sd                            block          
##    student_t(3, 0, 2.5)        sd                  Intercept block          
##    student_t(3, 0, 2.5)        sd                              vid          
##    student_t(3, 0, 2.5)        sd                  Intercept   vid          
##    student_t(3, 0, 2.5)     sigma                                           
##  nlpar bound       source
##                      user
##              (vectorized)
##              (vectorized)
##              (vectorized)
##              (vectorized)
##              (vectorized)
##              (vectorized)
##              (vectorized)
##              (vectorized)
##              (vectorized)
##              (vectorized)
##                   default
##                   default
##              (vectorized)
##              (vectorized)
##              (vectorized)
##              (vectorized)
##                   default
check_prior(consumption_assets_normal_bayesmodel)
##                       Parameter Prior_Quality
## 1                   b_Intercept   informative
## 2              b_cost_deviation uninformative
## 3                   b_treat_any   informative
## 4                    b_treat_GK   informative
## 5  b_assetsconsumption_asinh_R1 uninformative
## 6            b_Lhh_wealth_asinh uninformative
## 7        b_Lvill_eligible_ratio uninformative
## 8                  b_Lroomsnumb   informative
## 9        b_Ldurablesexpenditure uninformative
## 10    b_Ldietarydi_x_Lassetscon uninformative
## 11    b_Lproductiv_x_Lassetscon uninformative

Diagnostics

# trace diagnostic plot
mcmc_trace(consumption_assets_normal_bayesmodel, n_warmup = 0,
           pars = c("b_Intercept", "b_cost_deviation", "b_treat_any", "b_treat_GK",
                    "b_assetsconsumption_asinh_R1", "b_Lhh_wealth_asinh", "b_Lvill_eligible_ratio", 
                    "b_Lroomsnumb", "b_Ldurablesexpenditure",
                    "b_Ldietarydi_x_Lassetscon", "b_Lproductiv_x_Lassetscon",
                    "sd_block__Intercept", "sd_vid__Intercept", "sigma"))

ggsave("table_3_diagnostics\\consumption_assets_normal_race.png", plot = last_plot(), width = 12, height = 5)

#density diagnostic plots
mcmc_dens(consumption_assets_normal_bayesmodel,
           pars = c("b_Intercept", "b_cost_deviation", "b_treat_any", "b_treat_GK",
                    "b_assetsconsumption_asinh_R1", "b_Lhh_wealth_asinh", "b_Lvill_eligible_ratio", 
                    "b_Lroomsnumb", "b_Ldurablesexpenditure",
                    "b_Ldietarydi_x_Lassetscon", "b_Lproductiv_x_Lassetscon",
                    "sd_block__Intercept", "sd_vid__Intercept", "sigma"))

ggsave("table_3_diagnostics\\consumption_assets_normal_dens.png", plot = last_plot(), width = 12, height = 5)

mcmc_dens_overlay(consumption_assets_normal_bayesmodel,
           pars = c("b_Intercept", "b_cost_deviation", "b_treat_any", "b_treat_GK",
                    "b_assetsconsumption_asinh_R1", "b_Lhh_wealth_asinh", "b_Lvill_eligible_ratio", 
                    "b_Lroomsnumb", "b_Ldurablesexpenditure",
                    "b_Ldietarydi_x_Lassetscon", "b_Lproductiv_x_Lassetscon",
                    "sd_block__Intercept", "sd_vid__Intercept", "sigma"))

ggsave("table_3_diagnostics\\consumption_assets_normal_dens_overlay.png", plot = last_plot(), width = 12, height = 5)

#acf (auto-correlation) diagnostic plot
mcmc_acf(consumption_assets_normal_bayesmodel,
           pars = c("b_Intercept", "b_cost_deviation", "b_treat_any", "b_treat_GK",
                    "b_assetsconsumption_asinh_R1", "b_Lhh_wealth_asinh", "b_Lvill_eligible_ratio", 
                    "b_Lroomsnumb", "b_Ldurablesexpenditure",
                    "b_Ldietarydi_x_Lassetscon", "b_Lproductiv_x_Lassetscon",
                    "sd_block__Intercept", "sd_vid__Intercept", "sigma"))

ggsave("table_3_diagnostics\\consumption_assets_normal_acf.png", plot = last_plot(), width = 12, height = 5)

posterior predictive checks

pp_check(consumption_assets_normal_bayesmodel, ndraws = 100)

pp_check(consumption_assets_normal_bayesmodel, ndraws = 10, type = 'error_scatter_avg', alpha = .1)

House Value - Secondary Outcome

Load the house value model variables

house value: This is the basic benchmarking model utilizing a prior normal distribution with mean = 0 and 6 standard deviations. This model had a number of errors and the posterior estimates are unreliable.

dwelling_cost_normal_bayesmodel <-
  brm(formula = selfcostdwell_asinh | weights(samp_wgt) ~
        cost_deviation + treat_any + treat_GK +
        selfcostdwell_asinh_R1 + Lhh_wealth_asinh + Lvill_eligible_ratio + Lroomsnumb + Ldurablesexpenditure +
        (1 | block) + (1 | vid),
     data = dwelling_cost_data,
     family = gaussian("identity"),
     set_prior("normal(0,6)", class = "b"),
     seed = 1272022,
     warmup = 1000,
     iter = 2000,
     thin = 1,
     control = list(adapt_delta = .95, max_treedepth = 10),
     #backend = "cmdstanr",
     cores = 4, #overrides default 1 core
     #threads = 3,need to get cmdstanr package working here
     save_pars = save_pars(all = TRUE), # potentially allows for more post-processing functionality
     file = "informed_prior_outcomes\\dwelling_cost_normal_bayes")

Model Summery

summary(dwelling_cost_normal_bayesmodel)
## Warning: Parts of the model have not converged (some Rhats are > 1.05). Be
## careful when analysing the results! We recommend running more iterations and/or
## setting stronger priors.
## Warning: There were 2023 divergent transitions after warmup. Increasing
## adapt_delta above 0.95 may help. See http://mc-stan.org/misc/
## warnings.html#divergent-transitions-after-warmup
##  Family: gaussian 
##   Links: mu = identity; sigma = identity 
## Formula: selfcostdwell_asinh | weights(samp_wgt) ~ cost_deviation + treat_any + treat_GK + selfcostdwell_asinh_R1 + Lhh_wealth_asinh + Lvill_eligible_ratio + Lroomsnumb + Ldurablesexpenditure + (1 | block) + (1 | vid) 
##    Data: dwelling_cost_data (Number of observations: 1654) 
##   Draws: 4 chains, each with iter = 2000; warmup = 1000; thin = 1;
##          total post-warmup draws = 4000
## 
## Group-Level Effects: 
## ~block (Number of levels: 22) 
##               Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sd(Intercept)     0.83      0.66     0.20     1.93 3.32        4       11
## 
## ~vid (Number of levels: 247) 
##               Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sd(Intercept)     0.29      0.15     0.16     0.70 4.15        4       11
## 
## Population-Level Effects: 
##                        Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS
## Intercept                 12.56      8.41    -1.15    26.73 3.35        4
## cost_deviation             0.00      0.00    -0.00     0.00 2.59        5
## treat_any                  0.13      1.02    -1.45     1.51 2.44        5
## treat_GK                   1.51      0.27     0.94     1.79 3.17        4
## selfcostdwell_asinh_R1    -0.20      0.86    -1.33     0.80 2.56        5
## Lhh_wealth_asinh           0.07      0.06    -0.01     0.19 2.51        5
## Lvill_eligible_ratio       0.69      0.53    -0.19     1.28 3.42        4
## Lroomsnumb                 0.06      0.98    -1.43     1.31 3.07        5
## Ldurablesexpenditure       0.00      0.00     0.00     0.00 1.85        6
##                        Tail_ESS
## Intercept                    11
## cost_deviation               30
## treat_any                    22
## treat_GK                     11
## selfcostdwell_asinh_R1       31
## Lhh_wealth_asinh             22
## Lvill_eligible_ratio         11
## Lroomsnumb                   11
## Ldurablesexpenditure         21
## 
## Family Specific Parameters: 
##       Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sigma     3.06      1.95     0.89     8.63 2.91        5       11
## 
## Draws were sampled using sampling(NUTS). For each parameter, Bulk_ESS
## and Tail_ESS are effective sample size measures, and Rhat is the potential
## scale reduction factor on split chains (at convergence, Rhat = 1).

Prior summery - how informative are priors

prior_summary(dwelling_cost_normal_bayesmodel)
##                    prior     class                   coef group resp dpar nlpar
##              normal(0,6)         b                                             
##              normal(0,6)         b         cost_deviation                      
##              normal(0,6)         b   Ldurablesexpenditure                      
##              normal(0,6)         b       Lhh_wealth_asinh                      
##              normal(0,6)         b             Lroomsnumb                      
##              normal(0,6)         b   Lvill_eligible_ratio                      
##              normal(0,6)         b selfcostdwell_asinh_R1                      
##              normal(0,6)         b              treat_any                      
##              normal(0,6)         b               treat_GK                      
##  student_t(3, 13.8, 2.5) Intercept                                             
##     student_t(3, 0, 2.5)        sd                                             
##     student_t(3, 0, 2.5)        sd                        block                
##     student_t(3, 0, 2.5)        sd              Intercept block                
##     student_t(3, 0, 2.5)        sd                          vid                
##     student_t(3, 0, 2.5)        sd              Intercept   vid                
##     student_t(3, 0, 2.5)     sigma                                             
##  bound       source
##                user
##        (vectorized)
##        (vectorized)
##        (vectorized)
##        (vectorized)
##        (vectorized)
##        (vectorized)
##        (vectorized)
##        (vectorized)
##             default
##             default
##        (vectorized)
##        (vectorized)
##        (vectorized)
##        (vectorized)
##             default
check_prior(dwelling_cost_normal_bayesmodel)
##                  Parameter Prior_Quality
## 1              b_Intercept   informative
## 2         b_cost_deviation uninformative
## 3              b_treat_any   informative
## 4               b_treat_GK uninformative
## 5 b_selfcostdwell_asinh_R1   informative
## 6       b_Lhh_wealth_asinh uninformative
## 7   b_Lvill_eligible_ratio uninformative
## 8             b_Lroomsnumb   informative
## 9   b_Ldurablesexpenditure uninformative

Diagnostics

# trace diagnostic plot
mcmc_trace(dwelling_cost_normal_bayesmodel, n_warmup = 0,
           pars = c("b_Intercept", "b_cost_deviation", "b_treat_any", "b_treat_GK",
                    "b_selfcostdwell_asinh_R1", "b_Lhh_wealth_asinh", "b_Lvill_eligible_ratio", 
                    "b_Lroomsnumb", "b_Ldurablesexpenditure",
                    "sd_block__Intercept", "sd_vid__Intercept", "sigma"))

ggsave("table_3_diagnostics\\dwelling_cost_normal_trace.png", plot = last_plot(), width = 12, height = 5)

#density diagnostic plots
mcmc_dens(dwelling_cost_normal_bayesmodel,
           pars = c("b_Intercept", "b_cost_deviation", "b_treat_any", "b_treat_GK",
                    "b_selfcostdwell_asinh_R1", "b_Lhh_wealth_asinh", "b_Lvill_eligible_ratio", 
                    "b_Lroomsnumb", "b_Ldurablesexpenditure",
                    "sd_block__Intercept", "sd_vid__Intercept", "sigma"))

ggsave("table_3_diagnostics\\dwelling_cost_normal_dens.png", plot = last_plot(), width = 12, height = 5)

mcmc_dens_overlay(dwelling_cost_normal_bayesmodel,
           pars = c("b_Intercept", "b_cost_deviation", "b_treat_any", "b_treat_GK",
                    "b_selfcostdwell_asinh_R1", "b_Lhh_wealth_asinh", "b_Lvill_eligible_ratio", 
                    "b_Lroomsnumb", "b_Ldurablesexpenditure",
                    "sd_block__Intercept", "sd_vid__Intercept", "sigma"))

ggsave("table_3_diagnostics\\dwelling_cost_normal_dens_overlay.png", plot = last_plot(), width = 12, height = 5)

#acf (auto-correlation) diagnostic plot
mcmc_acf(dwelling_cost_normal_bayesmodel,
           pars = c("b_Intercept", "b_cost_deviation", "b_treat_any", "b_treat_GK",
                    "b_selfcostdwell_asinh_R1", "b_Lhh_wealth_asinh", "b_Lvill_eligible_ratio", 
                    "b_Lroomsnumb", "b_Ldurablesexpenditure",
                    "sd_block__Intercept", "sd_vid__Intercept", "sigma"))

ggsave("table_3_diagnostics\\dwelling_cost_normal_acf.png", plot = last_plot(), width = 12, height = 5)

posterior predictive checks

pp_check(dwelling_cost_normal_bayesmodel, ndraws = 100)

pp_check(dwelling_cost_normal_bayesmodel, ndraws = 10, type = 'error_scatter_avg', alpha = .1)

house quality - Secondary Outcome

Load the house quality model variables

house quality: This is the basic benchmarking model utilizing the default, uninformed priors

housing_quality_normal_bayesmodel <-
  brm(formula = housing_quality | weights(samp_wgt) ~
        cost_deviation + treat_any + treat_GK +
        housing_quality_R1 + Lhh_wealth_asinh + Lvill_eligible_ratio + Lroomsnumb +
        (1 | block) + (1 | vid),
     data = housing_quality_data,
     family = gaussian("identity"),
     set_prior("normal(0,6)", class = "b"),
     seed = 1272022,
     warmup = 1000,
     iter = 2000,
     thin = 1,
     control = list(adapt_delta = .95, max_treedepth = 10),
     #backend = "cmdstanr",
     cores = 4, #overrides default 1 core
     #threads = 3,need to get cmdstanr package working here
     save_pars = save_pars(all = TRUE), # potentially allows for more post-processing functionality
     file = "informed_prior_outcomes\\housing_quality_normal_bayes")

Model Summery

summary(housing_quality_normal_bayesmodel)
##  Family: gaussian 
##   Links: mu = identity; sigma = identity 
## Formula: housing_quality | weights(samp_wgt) ~ cost_deviation + treat_any + treat_GK + housing_quality_R1 + Lhh_wealth_asinh + Lvill_eligible_ratio + Lroomsnumb + (1 | block) + (1 | vid) 
##    Data: housing_quality_data (Number of observations: 1751) 
##   Draws: 4 chains, each with iter = 2000; warmup = 1000; thin = 1;
##          total post-warmup draws = 4000
## 
## Group-Level Effects: 
## ~block (Number of levels: 22) 
##               Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sd(Intercept)     0.11      0.07     0.01     0.27 1.01      635     1399
## 
## ~vid (Number of levels: 248) 
##               Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sd(Intercept)     0.59      0.05     0.50     0.68 1.00     2097     2836
## 
## Population-Level Effects: 
##                      Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS
## Intercept               -2.62      0.21    -3.02    -2.21 1.00     4181
## cost_deviation           0.00      0.00    -0.00     0.00 1.00     2716
## treat_any               -0.12      0.13    -0.38     0.14 1.00     2334
## treat_GK                -0.08      0.13    -0.34     0.18 1.00     2093
## housing_quality_R1       0.01      0.02    -0.03     0.05 1.00     5163
## Lhh_wealth_asinh         0.04      0.01     0.02     0.05 1.00     5704
## Lvill_eligible_ratio     0.14      0.43    -0.69     0.99 1.00     1752
## Lroomsnumb               0.52      0.04     0.45     0.59 1.00     4986
##                      Tail_ESS
## Intercept                3366
## cost_deviation           2751
## treat_any                2695
## treat_GK                 2609
## housing_quality_R1       3668
## Lhh_wealth_asinh         2826
## Lvill_eligible_ratio     2517
## Lroomsnumb               2988
## 
## Family Specific Parameters: 
##       Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sigma     1.92      0.02     1.88     1.96 1.00     6705     3269
## 
## Draws were sampled using sampling(NUTS). For each parameter, Bulk_ESS
## and Tail_ESS are effective sample size measures, and Rhat is the potential
## scale reduction factor on split chains (at convergence, Rhat = 1).

Prior summery - how informative are priors

prior_summary(housing_quality_normal_bayesmodel)
##                    prior     class                 coef group resp dpar nlpar
##              normal(0,6)         b                                           
##              normal(0,6)         b       cost_deviation                      
##              normal(0,6)         b   housing_quality_R1                      
##              normal(0,6)         b     Lhh_wealth_asinh                      
##              normal(0,6)         b           Lroomsnumb                      
##              normal(0,6)         b Lvill_eligible_ratio                      
##              normal(0,6)         b            treat_any                      
##              normal(0,6)         b             treat_GK                      
##  student_t(3, -0.1, 2.5) Intercept                                           
##     student_t(3, 0, 2.5)        sd                                           
##     student_t(3, 0, 2.5)        sd                      block                
##     student_t(3, 0, 2.5)        sd            Intercept block                
##     student_t(3, 0, 2.5)        sd                        vid                
##     student_t(3, 0, 2.5)        sd            Intercept   vid                
##     student_t(3, 0, 2.5)     sigma                                           
##  bound       source
##                user
##        (vectorized)
##        (vectorized)
##        (vectorized)
##        (vectorized)
##        (vectorized)
##        (vectorized)
##        (vectorized)
##             default
##             default
##        (vectorized)
##        (vectorized)
##        (vectorized)
##        (vectorized)
##             default
check_prior(housing_quality_normal_bayesmodel)
##                Parameter Prior_Quality
## 1            b_Intercept   informative
## 2       b_cost_deviation uninformative
## 3            b_treat_any uninformative
## 4             b_treat_GK uninformative
## 5   b_housing_quality_R1 uninformative
## 6     b_Lhh_wealth_asinh uninformative
## 7 b_Lvill_eligible_ratio uninformative
## 8           b_Lroomsnumb uninformative

Diagnostics

# trace diagnostic plot
mcmc_trace(housing_quality_normal_bayesmodel, n_warmup = 0,
           pars = c("b_Intercept", "b_cost_deviation", "b_treat_any", "b_treat_GK",
                    "b_housing_quality_R1", "b_Lhh_wealth_asinh", "b_Lvill_eligible_ratio", "b_Lroomsnumb",
                    "sd_block__Intercept", "sd_vid__Intercept", "sigma"))

ggsave("table_3_diagnostics\\housing_quality_normal_trace.png", plot = last_plot(), width = 12, height = 5)

#density diagnostic plots
mcmc_dens(housing_quality_normal_bayesmodel,
           pars = c("b_Intercept", "b_cost_deviation", "b_treat_any", "b_treat_GK",
                    "b_housing_quality_R1", "b_Lhh_wealth_asinh", "b_Lvill_eligible_ratio", "b_Lroomsnumb",
                    "sd_block__Intercept", "sd_vid__Intercept", "sigma"))

ggsave("table_3_diagnostics\\housing_quality_normal_dens.png", plot = last_plot(), width = 12, height = 5)

mcmc_dens_overlay(housing_quality_normal_bayesmodel,
           pars = c("b_Intercept", "b_cost_deviation", "b_treat_any", "b_treat_GK",
                    "b_housing_quality_R1", "b_Lhh_wealth_asinh", "b_Lvill_eligible_ratio", "b_Lroomsnumb",
                    "sd_block__Intercept", "sd_vid__Intercept", "sigma"))

ggsave("table_3_diagnostics\\housing_quality_normal_dens_overlay.png", plot = last_plot(), width = 12, height = 5)

#acf (auto-correlation) diagnostic plot
mcmc_acf(housing_quality_normal_bayesmodel,
           pars = c("b_Intercept", "b_cost_deviation", "b_treat_any", "b_treat_GK",
                    "b_housing_quality_R1", "b_Lhh_wealth_asinh", "b_Lvill_eligible_ratio", "b_Lroomsnumb",
                    "sd_block__Intercept", "sd_vid__Intercept", "sigma"))

ggsave("table_3_diagnostics\\housing_quality_normal_acf.png", plot = last_plot(), width = 12, height = 5)

posterior predictive checks

pp_check(housing_quality_normal_bayesmodel, ndraws = 100)

pp_check(housing_quality_normal_bayesmodel, ndraws = 10, type = 'error_scatter_avg', alpha = .1)